package tmhprediction.classification;

import java.io.File;
import java.io.IOException;
import java.util.ArrayList;
import java.util.Arrays;
import java.util.HashMap;
import java.util.LinkedList;
import java.util.List;
import java.util.Vector;

import org.apache.commons.io.FileUtils;

import util.FastaReader;

/**
 * Can achieve ~95% specificity and sensitivity, while processing an estimated 1000 sequences per second
 * @author jonas
 *
 */
public class SolDiscriminationBasedOnFreeInsertionEnergy {

	/**
	 * @param args
	 */
	
	static final int WINDOW_SIZE = 19;
	static final double CZERO = 2.70E-01;
//	static final double CONE = 9.29E+00;
//	static final double CTWO = 6.45E-01;
//	static final double CTHREE = 8.22E-03 ;
	
//	static final double CZERO = 0;
	static final double CONE = 0;
	static final double CTWO = 0;
	static final double CTHREE = 0;
	
	static double[][] energies = new double[20][WINDOW_SIZE];
	static HashMap<String, Integer> aaPos = new HashMap<String, Integer>();
	
	
	public static void main(String[] args)
	{
		String matrix = args[0];
//		String impf0 = args[1];
//		String impf1 = args[2];
//		String impf2 = args[3];
		
		String base = args[1];
		
//		String solf0 = args[4];
//		String solf1 = args[5];
//		String solf2 = args[6];
		
		readMatix(matrix);
		
		HashMap<String, String> impHash0 = FastaReader.readFastaFile(base+"/imp0.fasta");
		HashMap<String, String> impHash1 = FastaReader.readFastaFile(base+"/imp1.fasta");
		HashMap<String, String> impHash2 = FastaReader.readFastaFile(base+"/imp2.fasta");
		HashMap<String, String> solHash0 = FastaReader.readFastaFile(base+"/sol0.fasta");
		HashMap<String, String> solHash1 = FastaReader.readFastaFile(base+"/sol1.fasta");
		HashMap<String, String> solHash2 = FastaReader.readFastaFile(base+"/sol2.fasta");
		
		//Optimization
//		System.out.println("First");
//		evaluateThisPredictor(getValues(impHash1), getValues(solHash1));
//		System.out.println("Second");
//		evaluateThisPredictor(getValues(impHash2), getValues(solHash2));
//		System.out.println("Third");
//		evaluateThisPredictor(getValues(impHash0), getValues(solHash0));
		
		
//		//Test based on the results from above, manually evaluated
//		System.out.println("First");
//		evaluateThisPredictorGivenCutoffs(getValues(impHash2), getValues(solHash2),26 ,3.75);
//		System.out.println("Second");
//		evaluateThisPredictorGivenCutoffs(getValues(impHash0), getValues(solHash0),18,1.5);
//		System.out.println("Third");
//		evaluateThisPredictorGivenCutoffs(getValues(impHash1), getValues(solHash1),26,3.0);		
	
		
		
		//Optimization
//		System.out.println("First");
//		HashMap<String, String> impHash01 = impHash0;
//		impHash01.putAll(impHash1);
//		HashMap<String, String> solHash01 = solHash0;
//		solHash01.putAll(solHash1);
//		
//		evaluateThisPredictor(getValues(impHash01), getValues(solHash01));
//		
//		
//		HashMap<String, String> impHash12 = impHash1;
//		impHash12.putAll(impHash2);
//		HashMap<String, String> solHash12 = solHash1;
//		solHash12.putAll(solHash2);
//		System.out.println("Second");
//		evaluateThisPredictor(getValues(impHash12), getValues(solHash12));
//		
//		
//		-
//		System.out.println("Third");
//		HashMap<String, String> impHash20 = impHash2;
//		impHash20.putAll(impHash0);
//		HashMap<String, String> solHash20 = solHash2;
//		solHash20.putAll(solHash0);
//		evaluateThisPredictor(getValues(impHash20), getValues(solHash20));
		
		
		
		//Eval2
		evaluateThisPredictorGivenCutoffs(getValues(impHash2), getValues(solHash2),28,4.0);
		System.out.println("Second");
		evaluateThisPredictorGivenCutoffs(getValues(impHash0), getValues(solHash0),25,4.0);
		System.out.println("Third");
		evaluateThisPredictorGivenCutoffs(getValues(impHash1), getValues(solHash1),26,4.0);		
		
		
		
//		String sequence = "MPLRKMKIPFLLLFFLWEAESHAASRPNIILVMADDLGIGDPGCYGNKTIRTPNIDRLASGGVKLTQHLAASPLCTPSRAAFMTGRYPVRSGMASWSRTGVFLFTASSGGLPTDEITFAKLLKDQGYSTALIGKWHLGMSCHSKTDFCHHPLHHGFNYFYGISLTNLRDCKPGEGSVFTTGFKRLVFLPLQIVGVTLLTLAALNCLGLLHVPLGVFFSLLFLAALILTLFLGFLHYFRPLNCFMMRNYEIIQQPMSYDNLTQRLTVEAAQFIQRNTETPFLLVLSYLHVHTALFSSKDFAGKSQHGVYGDAVEEMDWSVGQILNLLDELRLANDTLIYFTSDQGAHVEEVSSKGEIHGGSNGIYKGGKANNWEGGIRVPGILRWPRVIQAGQKIDEPTSNMDIFPTVAKLAGAPLPEDRIIDGRDLMPLLEGKSQRSDHEFLFHYCNAYLNAVRWHPQNSTSIWKAFFFTPNFNPVGSNGCFATHVCFCFGSYVTHHDPPLLFDISKDPRERNPLTPASEPRFYEILKVMQEAADRHTQTLPEVPDQFSWNNFLWKPWLQLCCPSTGLSCQCDREKQDKRLSR";
		
//		String test = "MPLRKMKIPFLLLFFLWEA";
//		String[] bla = test.split("");
//		List<String> list = new ArrayList<String>(Arrays.asList(bla));
//		list.remove(0);
//		
//		System.out.println(calcEnergy(list));
		
		
//		double[] bla = calculateEnergiesForSequence(sequence);
//		for(double a :  bla)
//		{
//			System.out.println(a);
//		}
		
//		System.out.println(isIMP(sequence, 10));
		
	}
	
	private static List<String> getValues(HashMap<String, String> map)
	{
		List<String> result = new LinkedList<String>();
		for(String e : map.keySet())
		{
			result.add(map.get(e));//in theory redundancy
		}
		return result;
	}
	
	public static void evaluateThisPredictorGivenCutoffs(List<String> imps, List<String> sols, int i, double j)
	{
		Vector<Double> v = new Vector<Double>();
		v.add(Double.MIN_VALUE);
		int tp = 0, fp = 0 , fn = 0, tn = 0;
		
				for(String seq : imps)
				{
					if(isIMP(seq, i,j))
					{
						tp++;
//						countMissed++;
						
					}
					else
					{
						fn++;
					}
				}
				
				//all sol
				for(String seq : sols)
				{
					if(isIMP(seq, i,j))
					{
						fp++;
//						countOverpred++;
					}
					else
					{
						tn++;
					}
				}
//				System.out.println("Cutoff: " +i + " Cutoff Energy: " + j  + " Missed: " + countMissed + " Overpred: " + countOverpred);
				double sens = Double.valueOf(tp) / (tp + fn);
				double spec = Double.valueOf(tn) / (tn + fp);
				double combined = sens + spec;

					System.out.println("Cutoff: " +i + " Cutoff Energy: " + j  + " Sens: " + sens  + " Spec: " + spec + " combined: " + combined);
					System.out.println("TP: " + tp + " FP:" + fp + " FN:" + fn + " TN: " + tn);
	}
	
	public static void evaluateThisPredictor(List<String> imps, List<String> sols)
	{
		Vector<Double> v = new Vector<Double>();
		v.add(Double.MIN_VALUE);
		
		//try different cutoffs
		for(int i = 10; i < 40; i++)
		{
			for(double j = -3; j <= 5; j = j + 0.25 )
			{
//				int countOverpred = 0;
//				int countMissed = 0;
				
				int tp = 0, fp = 0 , fn = 0, tn = 0;
				
				//all imps first
				for(String seq : imps)
				{
					if(isIMP(seq, i,j))
					{
						tp++;
//						countMissed++;
						
					}
					else
					{
						fn++;
					}
				}
				
				//all sol
				for(String seq : sols)
				{
					if(isIMP(seq, i,j))
					{
						fp++;
//						countOverpred++;
					}
					else
					{
						tn++;
					}
				}
//				System.out.println("Cutoff: " +i + " Cutoff Energy: " + j  + " Missed: " + countMissed + " Overpred: " + countOverpred);
				double sens = Double.valueOf(tp) / (tp + fn);
				double spec = Double.valueOf(tn) / (tn + fp);
				double combined = sens + spec;
				if(v.lastElement() < combined)
				{
					v.add(combined);
					System.out.println("Cutoff: " +i + " Cutoff Energy: " + j  + " Sens: " + sens  + " Spec: " + spec + " combined: " + combined);
					System.out.println("TP: " + tp + " FP:" + fp + " FN:" + fn + " TN: " + tn);
				}
//				System.out.println("Cutoff: " +i + " Cutoff Energy: " + j  + " Sens: " + sens  + " Spec: " + spec + " combined: " + combined);
			}
		}
	}
	
	public static boolean isIMP(String seq, int cutoff, double  cutoffEnergy)
	{
		double[] arr = calculateEnergiesForSequence(seq);
		
		if(seq.length()<19)
		{
			throw new IllegalArgumentException();
		}
		int counter = 0;
		for(double d:  arr)
		{
			if(d < cutoffEnergy)
			{
				counter++;
			}
		}
//		System.out.println(counter);
		
		return (counter > cutoff) ? true : false;
	}
	
	public static double[] calculateEnergiesForSequence(String seq)
	{
		double[] energiesSeq = new double[seq.length()];
		//Divide the whole thing into substrings of length 19
//		if(seq.length()<19)
//		{
//			throw new IllegalArgumentException();
//		}
		
		for(int i = 0 ; i < seq.length()-18; i++)
		{
//			System.out.println(seq.substring(i, i+19));
			String[] split = seq.substring(i, i+19).split("");
			List<String> list = new ArrayList<String>(Arrays.asList(split));
			list.remove(0);
			
//			System.out.println(split[1]+" "+calcEnergy(list));
			energiesSeq[i] = calcEnergy(list);
		}
		return energiesSeq;
	}
	
	private static double getEnergyFromHash(String aa, int index)
	{
		double res = 0.0;
		if(aaPos.get(aa) != null)
		{
			res = energies[aaPos.get(aa)][index];
		}
		
		return res;
		
	}
	
	public static double calcEnergy(List<String> seq)
	{	
		//Calculate \sum_{i=1}^{l} \Delta G^{aa(i)}_{app}
		double sumOfPerPosEnergies = 0;
		for(int i = 0; i < seq.size(); i++)
		{
//			double energ = energies[aaPos.get(seq.get(i))][i];
			sumOfPerPosEnergies += getEnergyFromHash(seq.get(i), i);
		}
		
		//Calculate first factor under sqrt
		double firstFactorSum = 0;
		for(int i = 0; i < seq.size(); i++)
		{
//			double energ = energies[aaPos.get(seq.get(i))][i] + Math.sin(Math.toRadians(100.0*(i+1)));
			double energ = getEnergyFromHash(seq.get(i), i) + Math.sin(Math.toRadians(100.0*(i+1)));
			firstFactorSum += energ;
		}
		firstFactorSum = Math.pow(firstFactorSum,2.0);
		
		//Calculate second factor under sqrt
		double secondFactorSum = 0;
		for(int i = 0; i < seq.size(); i++)
		{
//			double energ = energies[aaPos.get(seq.get(i))][i] + Math.cos(Math.toRadians(100.0*(i+1)));
			double energ = getEnergyFromHash(seq.get(i), i) + Math.cos(Math.toRadians(100.0*(i+1)));
			secondFactorSum += energ;
		}
		secondFactorSum = Math.pow(secondFactorSum,2.0);
		
		//Complete root
		double root = Math.sqrt(firstFactorSum + secondFactorSum);
		
		double fullEnergy = sumOfPerPosEnergies + CZERO*(root) + CONE + CTWO*WINDOW_SIZE + CTHREE*(Math.pow(WINDOW_SIZE,2.0));
		
		return fullEnergy;
	}
	
	public static void readMatix(String file)
	{
		List<String> lines= null;
		try
		{
			lines = FileUtils.readLines(new File(file));
		} catch (IOException e)
		{
			e.printStackTrace();
		}
		
		int aacounter = 0;
		for (String line : lines)
		{
			String[] split = line.split(" ");
			String aa = split[0];
			
			aaPos.put(aa, aacounter);
			
			for(int i = 0; i < WINDOW_SIZE; i++)
			{
				energies[aacounter][i] = Double.valueOf(split[i+1]);
			}
			
			aacounter++;
		
		}
		
		
	}

}
